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The Formation of High Redshift Subminimeter Galaxies 
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ABSTRACT 

We describe a model for the formation of z ^ 2 Submillimeter Galaxies (SMGs) which 
simultaneously accounts for both average and bright SMGs while providing a reasonable 
match to their mean observed spectral energy distributions (SEDs). By coupling hydrody- 
namic simulations of galaxy mergers with the high resolution 3D polychromatic radiative 
transfer code SUNRISE, we find that a mass sequence of merger models which use observa- 
tional constraints as physical input naturally yield objects which exhibit black hole, bulge, and 
H2 gas masses similar to those observed in SMGs. The dominant drivers behind the 850 |J.m 
flux are the masses of the merging galaxies and the stellar birthcloud covering fraction. The 
most luminous (Ssso^lS mJy) sources are recovered by '^10^'^ Mq 1:1 major mergers with 
a birthcloud covering fraction close to unity, whereas more average SMGs (S85o~5-7 mJy) 
may be formed in lower mass halos (^5 x IO^^Mq ). These models demonstrate the need for 
high spatial resolution hydrodynamic and radiative transfer simulations in matching both the 
most luminous sources as well as the full SEDs of SMGs. While these models suggest a nat- 
ural formation mechanism for SMGs, they do not attempt to match cosmological statistics of 
galaxy populations; future efforts along this line will help ascertain the robustness of these 
models. 

Key words: cosmology:theory-galaxies:formation-galaxies:high-redshift- 
galaxies :interactions-galaxies : ISM-galaxies : starburst 



1 INTRODUCTION 

As a population, z ~ 2 Submillimeter Galaxies (SMGs) are the 
most lumino us, heavily star-forming galaxies in the Universe 
felain et al.ll2002h . The last decade of observations have provided 
significant constraints regarding their physical properties. Selected 
for their prodigious long wavelength emission ( >5 mJy at 850 
lJ.m), SMGs are thought to be starburst dominated ( Chapman et al. 
[2004; Alexander et al. 2008; Younger et al. 2008a), and have mas- 
sive~10^°--10"M f7) H2 gas reservoirs jTacconi et al j2()0^.l2008i : 
iGreve et alj20o3) . If stars form according to a Krou pa initial mass 
funct ion (IMF), as is supported by observations jTacconi et al.l 
I2OO8I) . then the infrared luminosities from SMGs translates to 
t remendous inferred ~700-1500 Mp yr~^ star formation rates 
Kovacs et al. I2OO6I; |PoDeetalJl20o3) riven by major mergers 
Tacconi et 'ZI I2OO8I) . MGs are extremely massive, with typical 
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~5x IO^^Mq bulges in place teorvs et al.ll2005l : [ Swinbank et al.l 
2004), and are a highly clustered populatio n of galaxies, re- 
siding in ^10^^-10" Mf7^ dark matter halos jBlain et al.1 12004| : 
ISwinbank etai]|2008h . 

Understanding SMGs from a theoretical standpoint has proven 
a major challenge for models. Recently, by assuming a flat stel- 
lar initial mass function (IMF) for starbursting galaxies (such that 
dn/d(ln m) = with x —Q as compared to 1.35 for Salpeter's 
semi-analytic models (SAMs) of galax y formation cou- 
pled with 2-D spectrophotometric simulations jSilva et ailll998h 
have successfully reproduced the observed SMG number counts at 
z ~ 2. Though these models have struggled to ma tch some details 
of SM Gs (e.g. the A'-band through mid-IR SEDs: ISwinbank et all 
|2008|) . it is clear from these models that, in principle, SMGs 



^ Though the iss ue of the IMF in z ~ 2 SMGs is controversial 
feaugh et al. 2005t IS winbank et al JI2008I: iTacconi et ai]|2008l:lDav3l2008l ; 
Ivan Dokkunil200^ . we do not attempt to investigate the role of the IMF in 
reproducing cosmological number counts in this paper; this will be investi- 
gated in due course. 
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broadly fit into our current understanding of hiigh-z galaxy forma- 
tion. 

However, while the SAMs provide insight into the cosmo- 
logical properties of SMGs, their usage of simplified analytic pre- 
scriptions for star formation and AGN activity mean that detailed 
information on the formation, evolution, and physical properties 
of individual SMGs is lacking. In this arena, high resolution hy- 
drodynamic simulations are better suited, yet thus far have strug- 
gled. While recent theoretical work utilizing gas rich galaxy merger 
simulations has shown good success in reproducing the observed 
characteristics of a number of galaxy populations - e.g. quasars 
jHopkins et al. l 2005, 200 ^, and refere nces therein), red galax- 
ies (Springel et al. 2005a; Hopkins et al. 2006b', '2008b'), elliptical 
galaxies dCox et al. 2006 ; Hopkins et aL ^ 08a.c, 2009a.h). and 
warm Ultraluminous Infrared Galaxies " jYounger et alj|2009al) - a 
detailed model for SMG formation and evolution using simulations 
has remained elusive. For example, attempts at forming SMGs by 
combining hydrodynamic and radiative transfer calculations of gas 
rich mergers have shown a peak flux of only ~4-5 mjy in the 
submillimeter for the most massive ( Mdm ~ 2x10^"^) and lumi- 
nous sources jChakrabarti et al.ll2008 >). This is only marginally de- 
tec tab le in current deep, wide cosmological surveys ( Coppin et al., 
l2006h . and has dark matter and stellar mas ses too large to be rep- 
resentative of the ave rage SMG population ( [Alexander et al.ll2008l : 
ISwinbank et"ai]|2008h . 

This Letter is the first in a series of papers attempting to under- 
stand the detailed properties of SMGs using simulations. Here, we 
present the first model for the formation of SMGs that reproduces 
both average and bright SMGs while plausibly matching their mean 
observed SEDs, and stellar bulge, H2, and black hole masses. We 
show that high resolution hydrodynamic and radiative transfer sim- 
ulations are necessary to capture both the requisite starbursts as 
well as stellar bulge growth and dust obscuration to simultane- 
ously power both the most luminous SMGs and match the observed 
SEDs. Throughout this paper, we assume Q.a = 0.7, Q.m = 0.3, 
and h = 0.7. 



2 NUMERICAL METHODS 

We follow the hydrodynamic evolution of galaxy mergers of vary- 
ing mass and mass ratios using the A'^-b ody plus smooth ed particle 
hydrodynamics (SPH) code GADGET3 l lSpringelll2005h . The syn- 
thetic photometric properties of the simulations are then calculated 
using the 3D adaptive m esh Monte Ca rlo polychromatic radiative 
transfer code, SUNRISE ( Jonssonl|2006l) . Our goal throughout is to 
utilize physical parameter input to the models directly constrained 
by observations of SMGs where possible, and from local ULIRGs 
or the Milky Way Galaxy when data from SMGs may not be avail- 
able. 



2.1 Hydrodynamics 

GADGET3 Utilizes an entr opy-conserving formalism for SPH 
jSpringel & Hernauisll2002|) which include prescriptions fo r radia- 
tive cooling of the gas jKatz etaLlll996l ; iDave et al.lll999l) . and a 
sub-grid formalism for a multi-phase interstellar medium (ISM) 
in whi ch cold clouds are embedd ed in a hot, pressure confining 
phase jSpringel & Hemquistll2003h . The cold clouds are allowed to 
grow through radiative cooling of the hot phase gas, and similarly, 
supemovae may evaporate the cold gas into the hotter phase. Nu- 
merically, supernovae pressurization of the ISM is handled through 



an effective equation of state (E OS). Here, we employ a stiff EOS 
with 5EOS = 1 (see Fig 4. of ISpringeletai]|2005bh which al- 
lows for stable disk evolution. Star formation in the simulations 
proceeds following a generalized Kennicutt- Schmidt (KS) relation 
(e.g. SFR oc p^'^) with normalization set to match the local sur- 
face density relation. This assumption is bolstered by observations 
of S MGs which hint at a KS index near the locally observed value, 
1.5 dBouche et al.ll2007h . 

Black holes are included in the simulations as sink parti- 
cles that accrete surrounding material following an Eddington lim- 
ited Bondi-Hoyle-Lyttleton paramaterization (e.g. iB ondi & Hovi3 
1944). A subresolution model for feedback is incorporated as an 
isotropic thermal coupling between the active galactic nucleus 
(AGN) and the surrounding ISM. The e fficiency of this coup ling 
is tuned to match the local M-a relatio n jDi Matteo et al.ll2005l) . 

The galaxies are initialized with a lHernauistl ( ll99Cl) dark mat- 
ter pro file, and virial properti es scaled to be appropriate for redshift 
z ( [Robertson et al.l2006h . The galaxies are initialized with cir- 
cular velocities Kirc = 320 500 k m s~^, motivated by t he cir- 
cular velocities for SMGs measured bv [Tacconi etall ( [20081) . This 

results in halo masses of A^dm ~1 0^^ 10^^ Mq, consistent 

with SMG clustering measurements ( [Blain et al. 2004). The galax- 
ies are bulgeless, and the initial gas fraction is 80%. This results 
in galaxies with gas fractions < 40% at final coalescence (owing 
to gas consumption by star formation) , comparable to estimates 
of observed SMGs (^ Bouche et al.|[200'7 ). Our gas particle masses 
were M = 1.4x10*' h~^MQ , and softening lengths 100 /i^^pc. 
In this paper, we consider 12 merger models varying mass, merger 
orbit and mass ratio. We consider 1:1 mergers in ~10^^ 5x10^^ 
and 1O^^M0 halos, and 1:3 and 1:12 mergers with a A^dm ~10 
Mq primary galaxy. 



2.2 Radiative Transfer 

SUNRISE considers the propagation of radiation between UV and 
millim eter wavelengths through a dusty medium (see [jonssonl 
( [2OO6I) for greater detail, as well as [jonsson et al. I ( [2009[) for im- 
portant updates, some of which are summarized here). SUNRISE 
tracks model photon 'packets' originated by sources - AGN, stel- 
lar clusters, and the ISM itself - as they traverse the dusty ISM 
using a Monte Carlo methodology. As in galaxies in nature, the 
flux emergent from the model galaxy/galaxies is determined by the 
photons that are emitted from the stars or AGN in the camera di- 
rection and escape without interacting and also the photons that 
scatter or are re-emitted by dust into the camera direction. Here, 
we utilized 8 different cameras distributed uniformly around the 
model galaxy/galaxies. SUNRISE treats dust self-absorption and re- 
emission in calculating the dust temperature. 

New for the simulations here, the AGN emits an empirical 
template SED deriv ed from observations of unobscured quasars 
( [Hopkins et ai][2007[) . assuming a radiative efficiency of 10%. The 
normalization of this input spectrum is set by the total bolometric 
luminosity of the central black hole. The spectru m emitted from 
stella r clusters is calculated utilizing Starburst 99 ( [Leitherer et al.[ 
[ 19991) . We assume a Kroupa IMF for the stars, whi ch is in good 
agreement with th e estimates for z ~ 2 galaxies ( [Davd [2OO8I; 
[Tacconi etal][2008[ : [van DokkunJ [2008t) . For young stellar clus- 
ters, which are still located in their nascent birthclouds of molec- 
ular gas, the input spectrum includes the effects of the HII re- 
gion s and photodissocia tion regions (PDRs) surrounding the clus- 
ters ( [Groves et al][2008[) . The models describe the evolution of the 
HII regions and PDRs analytically and use the photoionization 
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Figure 1. Evolution of 850 \im flux from model 1:1 major merger in a 
~10^^ Mq halo. Two different birthcloud covering fractions (fpdr) ^re 
shown (black solid line and red dashed line). The blue crosses represent the 
SFR with units on the right axis. The black star and dot near the top axis 
denote the peak of the starburst and black hole accretion rate, respectively. 
During the peak starburst, the opacity is dominated by stellar birthclouds, 
and the galaxy may be seen as an SMG. The most massive mergers with 
birthcloud covering fractions of unity may represent the most luminous ob- 
served SMGs. 



code MAPPINGSIII jGroves et aDl2004l) to calculate the SED that 
emerges from the cluster after it has been reprocessed by the HII 
regions and PDRs. The HII region absorbs almost all of the ion- 
izing radiation and is responsible for the hydrogen emission lines 
and hot dust emission. The PDR absorbs a significant fraction of 
the UV stellar light and reprocesses it to FIR emission via cooler 
dust. 

A fraction (fpoR) of clusters with ages < 10 Myr are mod- 
eled as obscured by their nascent birth clouds. As this fraction is 
increased, more UV and optical flux is reprocessed to longer wave- 
lengths. The covering fraction affects the maximum sub-mm lu- 
minosity during the peak of the starburst. This is equivalent to a 
cloud clearing time scale. While fpdr ~ 0.2 (which corresponds to 
a clearing timesca le ~2 Myr) m ay be a reasonable approximation 
for field galaxies ^Groves et al . 2008), in the case of massive, gas 
rich galaxy mergers, an assumption of fpDR ~ 1 (icicar >> 10 
Myr) may be more appropriate. Indeed, observations of molecular 
gas in local ULIRGs have found that the nuclear regions of mergers 
are often best characterized by a uniform molecular mediu m, rather 
than discrete clouds (Downes & Solomon 1998; Sakam oto et al.l 
1 1999) . While this is not the same as a birthcloud covering frac- 
tion of unity, tests done for this work have shown that they result 
in similar FIR/submm fluxe|f|. We note that the simulations do not 
account for the obscuration of stellar light by GMCs outside of the 
birthcloud. 

We utilize 10^ photon packets per iteration in the radiative 
equilibrium calculations, and the analysis is done for rest frame 
wavelengths 0.1 to 1000 /im. The stellar particles initialized with 
the simulation are assumed to have formed at a constant rate over 
~250 Myr. Metal enrichment is tracked from supemovae using an 

^ We have run experiments assuming fpd,. = 0, and that each giid cell is 
uniformly filled with dust with mass set by the ULIRG dust to gas ratio. In 
this limiting case, we find fluxes quite similar to those in simulations where 
the PDR covering fraction is unity. 
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Figure 2. Total halo mass versus peak Sgso. The full range of observed 850 
\xm fluxes from SMGs may be understood from a mass sequence merger 
models, varying only total mass. More massive mergers induce stronger 
bursts which fuel more luminous SMGs. The red star shows an isolated disk. 
Some merging activity may be necessaiy to induce the necessary SFRs to 
form SMGs. Even massive (~10^^Mq) isolated disks seem unable to reach 
the nominal Sgso > 5 mJy detection threshold, though may be routinely 
detectable via future-generation bolometers. 



instantaneous recycling approximation. The gas and stars initial- 
ized with the simulation are assigned a metallicity according to 
a closed-box model such that Z = (— y ln[/gas]) where Z is 
the metallicity, y the yield=0.02, and /gas is the initial gas frac- 
tion (though note the fluxes during the SMG phase of the model 
galaxies is not very sensitive to these assumptions). Without knowl- 
edge of the du st properties of z ~ 2 SMGs, we conservatively as- 
sume a 3. l_60 lDraine & Lil (l2007h model. We additionally assume 
a constant dust to gas ratio comparable to that of ULIRGs (1/50; 
Wilson et al. 2008) which tentative evidence suggests is appropri- 
ate for SMGs i Kovacs et al.ll2006l;lGreve et al.ll20()5l;lTacconi et all 
i2006i) . We have addition ally run mod els assuming a Milky Way 
dust to metals ratio of 0.4 ( lDwekll993) . and found similar results to 
within 10%. We note that comparisons to templates which include 
stochastically heated grains (Draine & Li 2007) suggest a potentia l 
uncertainly in the FIR SED up to a factor of 2 djonsson et^alj2009l) . 



3 RESULTS 

3.1 Formation of a Submillimeter Galaxy 

Galaxy mergers residing in halos with masses comparable to those 
observed naturally produce submillimeter galaxies with a range of 
850 \im fluxes during their peak starburst phase. Mergers with a 
total halo mass above '^5xlO^^M0 will produce galaxies above 
the nominal SMG detection threshold (Ssso ^ 5 mJy) whereas sig- 
nificantly lower mass mergers may not. The peak observed submil- 
limeter flux is dependent on the total dust obscuration for a given 
merger model. 

To see this, in Figure[T] we begin by showing the evolution of 
the 850 [i.m flux (at z = 2, the r ough mean redshift for radio-selected 
SMGs; IChapman et al.ll2003bh for a Mdm « 10^^ Mq 1:1 merger 
for two birthcloud covering fractions. Additionally, we overplot the 
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Table 1. 



Total Mdm 

Mq 


Mass Ratio 


Mq 


*%2 

Mq 


MbhI 

Mq 


Peak Ssso 
mjy 


2x10^^ 


1:1 


8x10" 


1.5x10" 


3.3x10** 


14.1 


1.1x1013 


1:3 


2.8x10" 


3x101" 




6.6 


7x10^2 


1:1 


2x10" 


1.8x10"' 


0.8x10** 


5.6 


9.8x1012 


1:12 


2.9x10" 


2x10"' 




4.8 



°' Owing to long dynamical friction time scales, minor mergers tend to be relatively gas 
poor upon final coalescence, and show their largest burst of 850 ^m flux on first passage. 
While BH growth is generally self-regulated, and the final coalescence BH masses are 
robust for a large range of initial BH see d masses, Mbh for fi rst passage starbursts can 
vary wildly based on initial seed masses iYounger et alj2008bl) . As such, the BH masses 
for first-passage SMGs in the simulations are relatively uncertain. 



SFR in light blue symbols. During the first passage and in-spiral 
phase of a galaxy merger, the galaxies are forming stars (individu- 
ally) at ~ 100-200 MQyr^^. During this time, obscuration by the 
diffuse ISM dominates and produces a modest ~5 mJy of emission 
seen during the in-spiral phase, regardless of the GMC covering 
fraction. 

As the galaxies reach final coalescence (T~0.65 Gyr), tidal 
torquing of the gas fuels a ~2000 MQyr^^ nuclear starburst. At 
this time, the FIR/submillimeter flux is dominated by cold dust re- 
processing of UV photons in young stellar clusters by their birth- 
clouds. Hence, the submillimeter flux rises rapidly with this burst 
in star formation rendering the galaxy detectable as an extremely 
luminous SMG for a relatively short (<50 Myr) lifetime. The peak 
flux of the SMG produced in the final burst is dependent on the 
birthcloud covering fraction. The merger model in Figure[T]with a 
covering fraction of 30% reaches a peak flux of ~7-8 mJy, whereas 
the model with covering fraction of unity produces an SMG com- 
parable to the most luminous observed sources. The obscuration of 
young stars is vital to creating the most luminous SMGs as these 
stars dominate the bolometric luminosity of the galaxy during coa- 
lescence. 

Because the peak flux received is primarily determined by the 
strength of the starburst (for a given birthcloud covering fraction), 
the submillimeter flux is dependent on the mass of the merger. Less 
massive galaxies will induce smaller tidal torques on the gas, fuel- 
ing lower SFRs; in these cases, the lightcurve shown for the mas- 
sive merger in Figure [T] will retain its shape, though simply lower 
in normalization. In Figure (2] we indicate this more explicitly by 
showing the predicted 850 \im flux against the galaxy (halo) mass 
of all 12 merger models varying primary halo mass, mass ratio (1:1, 
1:3 and 1:12) and orbit. Again, the 850 \im flux is modeled &t z = 
2. A single isolated disk is shown by the red star. Each is assumed 
to have a constant birthcloud covering fraction of unity. 

The 1:1 major mergers in the most massive (A/dm ~ 2 x 
IO^'^Mq ) halos produce SMGs comparable to t he most luminous 
observed jCoppin et al\ l2006l : IPope et alj |2006|) . Galaxies ~2^ 
times smaller may be representative of more average (Sgso ~ 5 
mJy) SMGs. Galaxies well below this mass limit are unlikely to 
ever produce SMGs based on the current fiducial criteria Sgso ^J5 
mJy, though may be detectable with sensitive new bolometer ar- 
rays. As is shown by the isolated galaxy (red star, Figure|2}, some 
merging seems to be necessary to fuel the starbursts that power 
SMGs. Isolated galaxies residing in even the most massive halos 



in our simulations (A/dm ~ 10^ Mq ) never appear to get above 
'-^2-3 mJy at 850 (xm. 



3.2 Physical Properties and SEDs of SMGs 

The SMGs produced by the galaxy mergers presented here have 
black hole, stellar, H2 , and dark matter masses (the latter by con- 
struction) comparable to those observed. In Table [S] we show 
these values for a typical orbit (6\,cf)\ = (30,60)°, 62,4)2 = 
(—30, 45)°) for the simulations which produced detectable SMGs. 
As shown in Cox et al. (2009, submitted) and [Narayanan et al.l 
I2OO8), binary galaxy mergers produce ~90% of their final stel- 
lar mass during the pre-coalescence phase. As such, the bulk of the 
~10"Mq bul ges are already in place during when the galaxy is 
visible as an SMG. Indeed, this is clear from a visual integration of 
the SFR in Figure[TJ an SFR of a fewx 10^ MQyr"^ for ~5 x 10* 
years results in a fewx 10^^ Mq bulge. 

The black hole growth follows a different story: the bulk of the 
growth occurs concomitant to the SMG phase/final coalescence. 
While the final black hole masses in the most massive/luminous 
systems may be similar to those of 2: ~ 2 quasars (~10^ Mq ), 

during the SMG phase they can be a factor of ~ 3 10 lower 

than these final values. This give s them masses of order ~10*M(^ , 
comparable to observed values ( [Alexander etai]|2008l and refer- 
ences therein). Finally, as was noted in §|2l by starting with pro- 
genitor galaxies with gas fractions fg = 0.8, the final gas fraction 
of during the SMG phase is <40% , comparable to estimates of ob- 
served SMGs teouche et al.l2007l) . While GADGET3 does not track 
the evolution of molecular gas, we derive the molecular gas frac- 
tion in post-processing based on the ambient pres sure on the cold, 
star-forming ISM (see iBlitz & RosolowskvluOOd ). From this, we 
find total H2 gas masses between 10^° - -lO" Mp (Tab le O. 
This is comparable to measurements by iGreve et alT il 20051) and 
iTacconi et al. Moreover, complementary radiative transfer 

simulations by|Naravanan et al. (2009) utilizing this methodology 
for assigning the H2 gas fraction have shown good correspon- 
dence between the simulated CO properties (e.g. excitations and 
line widths) of these model SMGs, and those in nature. 

Galaxy mergers in a mass range of Mdm ~5 x IO^^Mq - 
IO^'^Mq naturally produce SMGs with optical through mm-wave 
SEDs comparable to those observed. In Figure [3] we show the 
model SED for all snapshots in the merger simulations studied here 
which produce an SMG with Sgso > 5 rnjy, redshifted to z = 2. The 
blue shaded region shows the la dispersion amongst all snapshots 
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which satisfy this fiducial selection criterion while the black solid 
line shows the mean amo ngst snapshots. The m ean observed data 
points from the surveys of lChapman et al.' (2005. purple triangles), 
Hainline et al. (2009, red crosses), and Kovacs et al. (2006, green 
diamonds) are overlaid. Only observational data with spectroscopic 
redshifts are utilized. 

The SED encodes information regarding the young stars (at 
2: ~ 2, observed-frame A'-band), stellar mass (observed mid-IR), 
and cold dust obscuration levels (observed submillimeter). The 
SFR and stellar bulge growth are explicitly tracked in the progeni- 
tor galaxies and SMG by the high resolution SPH simulations while 
the dust obscuration is modeled globally by diffuse dust, and on 
local scales by the subgrid implementation of the MAPPINGSIII 
photoionization calculations. Because the SFRs in the simulated 
galaxies are comparable to those observed, and bulges of order '--^ 1~ 
8x10" Mq are in place during the SMG phase (see Table O, the 
observed if-band and mid-IR fluxes in the model SED match the 
observed data points closely. Sufficient cold dust obscuration pro- 
vided by the diffuse dust and birthclouds reprocesses the rest-frame 
UV emission into the cold dust tail providing a good match to the 
mean observed 350, 850, and 1 100 |J.m data points. The cases with 
lower birthcloud covering fraction (fpdr ~ 0.3) tend to overesti- 
mate the rest-frame UV and optical flux by a factor of ~2. 

While our simulations do not include radio emission, we can 
estimate the radio properties of our simulations by assuming that 
the f ar-IR/radio correlation does not evolve strongly with redshift 
(e.g. iKovacs et al]|2006l : I Younger et al.ll2009bl : ISaiina et alj|2008l 
Younger et al. 2009, in press). We estimate the FIR flux via the 
standard FIR = 1.26 x 10""(2.58S'6o+5'ioo) W m-^ where Seo 
and Sioo are the observed-frame fluxes at 60 and 100 [im, in Jy), 
and then relate the FIR flux to the radio (1.4 G Hz) flux through th e 
FIR-radio correlation (e.g. Equation 5 from iKovacs et al. 1 1200^ . 
We assume a range of (/-parameters, from 2.14 (consistent with ob- 
servations of SMGs; IKovacs et al.ll2006l), to the loc allv-observed 
value for IRAS BGS galaxies, q = 2.35 jSanders & M irabel 1996). 
We include bandwidth compression and a A-correction (assum- 
ing an SED slope of a — 0.7, typical for synchrotron power 
laws; ICondoiJ[l992h when redshifting the radio flux. We note 
that this method neglects any potential radio contribution from the 
AGN. When using the observed FIR-radio correlation for SMGs, 
we recover inferred 1.4 GHz flux densities (at z = 2.0) of 20- 
700 fiJy. The average value amongst our simulation sample for a 
g=2. 14(2.35) is 5i.4 ~ 125(95) /xJy, similar to values described 
in the compilation bv lChapman et al.l ( l2003al) . The 850/S'i.4 flux 
density ratios from our simulated SMGs is additional ly in good 
agree ment with measurements from recent surveys (e.g. IPope et aH 
I2OO6I) . For example, we find the mean 850/5*1.4 ratio for Ssso 
> 5mJy SMGs to be ~60 (with a la dispersion of 40) when utiliz- 
ing g=2.14. The mean(dispersion) increases to 100(65) when uti- 
lizing the local q. 

The SMGs formed in these simulations would generally have 
a high detection fraction of radio counterparts. In Figure|4] we plot 
the inferred 1.4 GHz flux density (derived utilizing q = 2.14) ver- 
sus Sgso at z=2.0 for all model galaxies in our simulations. We 
extend the nominal detection threshold to S85o> ImJy to serve 
as a prediction. SMGs (S850> 5 mjy) at z=2.0 will typically be 
detectable down to ~100 ^Jy at 1.4 GHz. While it is difficult to 
make quantitative comparisons to observed detection fractions, ow- 
ing to the fact that these simulations are not cosmological (and 
Ssso is roughly constant with increasing redshift while 5i.4 de- 
creases rapidly with increasing z), it is still feasible to make qual- 
itative comparisons to observations. Observed SMGs have a high 
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Figure 3. Mean SED of all snapshots from merger simulations satisfying 
the Sg5o>5 mJy selection criterion for SMGs (solid line). The dispersion 
amongst the snapshots is shaded in blue (the dispersion amongst viewing 
angles, however, is not shown). Th e simulated SED has b een redshifted to 
z = 2. Observed data points from IChaprnan et al] j2005l) . Hainline et al. 
(2009, submitted) and iKovacs et all j2006l) are overlaid (purple triangles, 
red crosses and green diamonds, respectively). These data all have spectro- 
scopic redshifts. The points represent the mean observed value from these 
surveys, though galaxies with 2 < 1 and upper limits have been discarded. 
The model SED from the tiducial merger provides a good match to the ob- 
served photometric data from the B-band to 1 mm, making this the first 
model SMG to reproduce the full optical through mm-wave SED. 



incidence of radio counterpa rts. Large radio ID fractions were re- 
ported bv llvison et al.l ( [2002h . who found that nearly all Ssso > 10 
mJy SMGs were detected above ~35 /xJy at 1.4 GHz, and approxi- 
mately half the SMGs above Ssso > 5 mJy were radio-detected. 
More recently. Wage et al. (2009) find that some ~80% of the 
SMGs in their sample contain radio counterparts above ~25/iJy. 
More generally, the dynamic range probed by the submillimeter 
flux densities and radio flux densities in Figure [4|is comparable 
to that seen in the S8so-5'i.4 relation compiled by Chapman et al.l 
( l2003al) . 

Many of the model SMGs in these simulations have multiple 
radio counterparts. In Figure[5] we show the optical (SDSS z, r, u) 
morphologies of our most luminous simulated SMG through its in- 
spiral with a contour denoting the location of the bulk of the ra- 
dio emission (contours enclose regions where Si. 4 > 0.1 ^Jy). 
As the galaxies in Figure |5] spiral in toward coalescence, they may 
be visible as multiple radio sources. There is a trend for the more 
massive/luminous SMGs to have multiple counterparts for a larger 
fraction of their lives. This owes to the fact that the most massive 
mergers (e.g. Mdm > 10" M0) can be SMGs during their in- 
spiral phase (e.g. Figure [T), whereas less massive galaxies (e.g. 
Mum ~ 5 X 10^^ Mq) generally only reach Sgso > 5 mJy when 
the progenitor galaxies have coalesced. To identify individual radio 
sources, we search for the peaks in the radio emission (as inferred 
from the FIR-radio correlation) in maps which view the galaxies 
at an arbitrary angle. We then assign the peaks as originating in 
separate galaxies when they originate at least 5 kpc apart. Utilizing 
this prescription, we find that the least massive/luminous SMGs in 
our simulations never have multiple radio counterparts. However, 
the SMGs formed in Mdm « 1(2) x 10^^ Mq halos have mul- 
tiple radio counterparts ~50(75)% of the time that they would be 
submillimeter-identified (Sg5o> 5 mJy). 
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Figure 5. Optical (Sloan 2, r, u) colors mapped onto R,G,B for the most massive/luminous SMG in our simulation sample. The thick green contours toward 
the galactic nuclei of each image enclose the region where the radio emission is above 5i.4ghz > 0.1 fiJy. Note, the contours are only at a single level. 
During the inspiral phase of this massive (Mjjm ~ 2 x 10^"^ 1^0) merger, the galaxy may be viewed as having multiple radio counterparts. 




100 

1.4 GHz Flux (nJy) 

Figure 4. Sgso versus 1.4 GHz inferred flux density at 2=2.0 for all models 
in Table |3] Given current flux limits, the SMGs simulated here will have a 
high incidence of radio counterparts. 



4 DISCUSSION AND CONCLUSION 

The models presented here have shown that merger-driven SMG 
formation simulations which utilize observational parameters as in- 
put produce both average (Sgso^S mjy) and bright (Ssso-^ 15 mjy) 
SMGs while reproducing observed physical characteristics, as well 
as the full optical-mm wave SED. This match between simulations 
and observations is only possible because of the high resolution ap- 
proach employed for both the hydrodynamic and subgrid physics 
included in the radiative transfer simulations. It is worth briefly 
comparing these results to previous efforts in modeling SMG for- 
mation to place our models in context. 

Models attempting to form SMGs bv lChakrabarti et al.l jZOOSi) 
utilized binary galaxy mergers similar to these. The peak submil- 
limeter flux from mergers comparable to the most massive observed 
(Mdm Ri2xlO^^MQ) in this work showed ~4-5 mJy at 850 nm; 
this is compared to ~15 mJy from comparable mass simulations at 
an identical merger orbit presented here. While the hydrodynamic 
evolution of the galaxy mergers is similar in both studies, the ra- 
diative transfer differs significantly. In particular, here, we employ 
a subgrid model for including the obscuration by the nascent birth- 
clouds of young stellar clusters. As Figure [T] shows, this obscura- 



tionby cold dust (which was not included in the lChakrabarti et^ 
(2008) simulations) is vital during the final burst to reprocess UV 
photons from young stars into FIR/submillimeter light[j 

The predicted SEDs from SMGs formed by semi-analytic ef- 
forts show underluminous Jf-band and mid-IR fluxes by up to an 
order of magnitude. This may owe to SFRs and stellar masses in 
these SMGs which are roughly an order of magn i tude l ower than 
those inferred by observations iBaugh et al Jl 20051 [2007h . While it 
is difficult to directly compare semi-analytic prescriptions with full 
A'^-body/SPH simulations, it is clear that varying treatments in the 
gas physics result in discrepant star formation histories, and con- 
sequent bulge masses and simulated SEDs. The plausible corre- 
spondence between the physical parameters and model SED from 
the models presented here and those observed owes to the detailed 
tracking of the starburst, star formation history, dust geometry and 
density distribution by these simulations. This suggests high reso- 
lution simulations are necessary to model the physical evolution of 
SMGs in detail. 

That said, while the simulations presented in this work ap- 
pear to provide a reasonable match to many observed parameters 
from SMGs, it is important to note that they, unlike the SAMs, are 
not cosmological in nature. For a galaxy formation model to be 
generic, it must concomitantly match the detailed physical and ob- 
servable properties of a given population, as well as cosmological 
statistics. In this sense, the simulations presented here are com- 
plementary to the SAMs of SMG formation. The former repro- 
duce detailed observables of SMGs in isolation via an ab initio 
model for SMG formation, though provide no information regard- 
ing cosmological statistics. The latter has some difficulty reproduc- 
ing e.g. the observed SED, though accurately matches the observed 
number counts (as well as present day B, K and 60M.m-ban d lu- 
minosity functions: iBaugh et al.|[2005l ; ISwinbank et alj|2008h . The 
ideal scenario would be a coupling of the methods by convolv- 
ing the lightcurves from high resolution hydrodynamic simula- 
tions of_2ndivMiM mergers to cosmological galaxy merger rates 
(e.g. iHopkins et aLll2009ch to investigate the statistical properties 
of SMGs. These efforts are currently are underway. 



^ IChakrabaili et alj )2008h incorporate a simpHfied, analytic model for the 
obscuring dust surrounding young stars. However, in contrast to our results, 
they find surprisingly Utile impact on the emergent IR SED, and even a 
moderate increase the rest-frame optical and near-UV luminosity. 
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